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Magnetic resonance imaging method 

[0001] The present invention relates to a magnetic resonance imaging method 
wherein undersampled magnetic resonance signals are acquired by a receiver 
antenna system having spatial sensitivity profiles according to the preamble portion 
of claim 1. The invention is further directed to a magnetic resonance apparatus and a 
computer program product for executing the magnetic resonance method. 
[0002] Such a magnetic resonance imaging method is usually indicated as a 
parallel imaging method and known from the paper by K. Pruessmann et al. in 
Magn. Reson. Med. 42 (1999), pages 952-962 and indicated as SENSE. Parallel 
imaging (PI) techniques exploit the spatial encoding effect inherent to 
inhomogeneous receiver sensitivity. Using arrays of multiple receiver coils, parallel 
acquisition permits reducing the /c-space sampling density, thus enabling significant 
scan time reduction. PI reconstruction according to the SENSE method is typically 
performed by image domain unfolding, which enables very fast processing in the 
common case of Cartesian /c-space sampling (Cartesian SENSE). The efficiency of 
this approach stems from the underlying assumption that the coil sensitivities do not 
change significantly across each voxel. While well fulfilled at high resolution, this 
assumption is gradually violated at low resolution, such as in spectroscopic imaging 
(SI). As a consequence, the spatial response function (SRF) in conventional SENSE 
spectroscopic imaging is prone to imperfections, which typically result in residual 
aliasing artifacts. 

[0003] To resolve the aliasing, Cartesian SENSE reconstructs each set of 
aliased voxels separately, thus decomposing the reconstruction into a series of 
smaller independent problems, which can be solved quickly. However, with this 
decomposition, the aliasing is only resolved exactly at the voxel centers. Residual 
aliasing artifacts occur when the sensitivity near a coil rises at a faster rate than the 
decay of the sine function. These artifacts are accentuated if /c-space is partially 
zero-filled. Thus, residual artifacts can arise in high-resolution MRI, even with 50% 
zerofilling in /c-space, which is often used in 3D applications. 

[0004] This problem was previously addressed by sensitivity extrapolation (see. 
U. Dydak et.al., Magn. Reson. Med. 2001; 46: p.71 3-722) and by the transition to 
modified, effective sensitivity values in a recursive scheme (see X. Zhao et. al., Proc. 
ISMRM 2002, Abstract No. 2393) 



2 

[0005] It is an object of the present invention to eliminate in the above 
mentioned parallel imaging method the residual aliasing artifacts arising from high 
resolution images reconstructed with parallel imaging techniques. It is a further 
object of the invention to provide a magnetic resonance imaging system and a 
computer program product adapted to such a method. 

[0006] This object is solved by the method as claimed in claim 1. The further 
objects are achieved by the magnetic resonance imaging system according to claim 
14 and the computer program product according to claim 15. 

[0007] The gist of the present invention is that residual aliasing artifacts can be 
eliminated by the minimum-norm least-squares reconstruction, which is consistent 
with the "strong" reconstruction as described in the above mentioned article of K.P. 
Prussmann et. al. 

[0008] Further advantages of the invention are disclosed in the dependent 
claims and in the following description in which an exemplified embodiment of the 
invention is described with respect to the accompanying drawings. It shows 
[0009] Fig. 1 an aliased image and the spatial response function of the 
marked pixel in the image, 

[0010] Fig. 2 a minimum-norm reconstruction at different iterations of the 
conjugate gradient algorithm, 

[0011] Fig. 3 the panels of a zoomed portion of the reconstructions, 
[0012] Fig. 4 SRF with minimum norm and conventional reconstruction, 
[0013] Fig. 5 phantom, conventional and minimum-norm reconstruction, and 
[0014] Fig. 6 diagrammatically a magnetic resonance imaging system 

General principles of the invention 

[0015] Consider a general MR imaging experiment during which complex 
data samples are acquired. Each measured signal value represents the spatial 
integral of an object-specific signal density function d(r) , modulated by an 
encoding function ^(r) : 
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where and denotes the complex noise component acquired 

along with the resonance signal. Noise in MRI is typically of thermal original, 
stemming from within the object as well as from the detector circuitry. 
[0016] The encoding functions can reflect any modulation that the object's 
magnetization undergoes in the course of the experiment. Typical means of spatial 
encoding are plane-wave modulation by linear gradient fields, sensitivity encoding by 
one or multiple receiver coils, or RF encoding by RF pulses. The encoding functions 
can also account for undesired effects, such as phase errors due to main field 
inhomogeneity. Incorporating only gradient and sensitivity encoding, for instance, the 
encoding functions read 

e»(r) = e"' V ^(>0 [2] 



where denotes the position in /c-space at which the value is taken and 

s ti( r ) denotes the complex spatial sensitivity of the receiver coil with which it is 
taken. 

[0017] Reconstructing an image from the acquired data amounts to calculating 
the values of the image components, which are typically called pixels. Let the index 

Jt count the pixel values p K in the desired resulting image, with N a denoting their 
total number. Typically, the calculation of the pixel values is assumed to be linear, 
reflecting the linearity of the integral in Eq. [1] with respect to d(r) . Then image 
reconstruction can be written as 

p = Fm [3] 
where the pixel values and the measured data were assembled in the vectors p 
and m , respectively, which are connected through the reconstruction matrix F with 
N x rows and columns. 

[0018] The problem of image reconstruction is thus recast as the problem of 
choosing the entries of the matrix F . The main goal in the choice of F is to ensure 
the proper depiction of the object, i.e. to make each pixel value reflect signal that 
originates from a corresponding small volume (voxel) in real space. In order to track 
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the spatial weighting of MR signal in the pixel value p n , one can insert Eq. [1] into 
Eq. [3] and extract the corresponding row: 

[0019] Here the expression in brackets reflects the spatial weighting of MR 
signal in the pixel value and is thus called the spatial response function (SRF) of the 
pixel: 

[5] 

[0020] Each pixel in the final image has an individual SRF. In order to make a 
pixel represent resonance signal from the corresponding voxel, its SRF should 
ideally be non-zero inside the voxel and zero outside the voxel. However, such ideal 
SRFs can usually not be accomplished, as they must be linear combinations of the 
encoding functions (see Eq. [5]). Therefore, SRFs in MRI are generally imperfect in 
that they exhibit smooth rather than sharp transitions at the voxel borders and are 
not exactly zero outside the voxel. In the resulting image, these defects are 
perceived as artifacts (e.g. blurring or aliasing), which limit the image's usefulness, 
e.g., its diagnostic value. 

[0021] In response to this problem, it is one object of the present invention to 
optimize the SRF, either globally or on a pixel by pixel basis. For this purpose, the 
SRF is preferably discretised by sampling the encoding functions along a finite set of 
positions in 3D space. Let this set be given by 

where P e fc~",} [6] 
[0022] Typically this set would form a regular Cartesian grid, yet it does not 
need to. The discretisation may as well be performed along any finite set of 
positions, e.g., along hexagonal or radial grids. Sampling the encoding functions at 

r p yields the encoding matrix 

which has rows and N p columns. Using this notion, the discretised SRF of the 
pixel Jt is given by the Jt -th row of the SRF matrix 

SRF discrete = FE [8] 
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[0023] More generally, each row of this matrix represents the discretised SRF 
of the corresponding voxel. Note that this representation can be taken to any level of 

accuracy by enhancing the number and density of the sampling positions r p . 

[0024] For SRF optimisation, the entries of F are chosen such that the 
individual SRFs approach given idealised target shapes. It depends on the 
application which SRF shape is considered ideal. In terms of overall depiction 
fidelity, Dirac distributions or boxcar functions may be regarded as optimal. For the 
purpose of suppressing far-range signal contamination, Gaussian shapes may be 
superior or potential prior knowledge about the resonance signal distribution in the 
object may favour other, tailored shapes. Any target SRF shape can be 
accommodated within the formalism explained here. Moreover, each pixel can be 
assigned an individual target SRF shape. This is done by a target SRF matrix T t 
whose 7r-th row is the discrete representation of the SRF considered ideal for the 
pixel with index tz\ 

T„, p =SRF? r * et (r p ) [9 ] 
[0025] The difference between the target SRFs and those actually achieved is 
hence reflected by ^RF discrete - T . For pixel-by-pixel SRF optimization, this 
difference is minimized row by row. In order to do so, the magnitude of each row of 
SRF d iscrete ~T \s typically expressed using a vector norm along with a general linear 
mapping preceding the norm operation. The SRF deviation for the pixel it is then 

L K -\((FE-T)A\\ [10] 

where A is a matrix representing a suitable linear mapping, the subscript Jt on the 
right side indicates the Jr-th row vector of the matrix product in brackets, and ||-|| 

denotes the chosen vector norm. In principle, any norm can be used, including the 
sum, maximum, and Euclidean norms. The Euclidean norm ||-|| is preferable as it 
leads to particularly feasible calculations. The SRF deviation is then quantified as 

A„-l((FE-T)Al\\ 2 -^((FE-T)e(FE-Tf)^ m 
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where 0 = AA H and the superscript H indicates the Hermitean adjoint, i.e. the 
complex conjugate transpose. 

[0026] Minimizing A^ independently for each musing common means of linear 

algebra yields the reconstruction matrix 

F = TdE H (EdE H Y [12] 

where the superscript + indicates the Moore-Penrose pseudoinverse. 
[0027] Reconstruction according to Eqs. [3, 12] yields optimal SRFs. One 
problem that may occur in applying this formula is ill-conditioning of the matrix F. Ill- 
conditioning leads to enhancement of the noise component in Eq. [4] (second term 
on right hand side of Eq. [4]), which frequently reaches unacceptable levels. 
Therefore, another object of this invention is to control noise enhancement in 
addition to SRF quality. According to another paper by KP Pruessmann et al. (Magn 
Reson Med 2001;46:638-651), the noise level in the pixel value x is given by 



* ^ V [13] 
where W denotes the noise covariance matrix of the acquired data. The properties 
of this matrix and the means of determining it are described in the above mentioned 
paper by KP Pruessmann et al. 

[0028] Using these notions, the task of ensuring proper image quality is 

restated as choosing F such that both and A^ are close to minimal for each 

pixel jr. This can, for instance, be accomplished by minimizing the weighted sum of 
squares of the two terms: 

min(a A 2 „ +fi A^ ) [14] 

where a and are scalar, real-valued, non-negative weighting coefficients. 

Performing this minimization using common means of linear algebra yields the 
reconstruction matrix 

F = fi T9E H (f$E0E H +aW) + [15] 

For a * 0 and /3 ^ 0 , an equivalent formula is given by 

F = a- l T(a- l E H W + E + fi' l e*y E H W + [16] 
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[0029] The equations [3, 15, 16] determine image reconstruction that controls 
SRF quality and noise simultaneously, individually for each pixel. To use these 
formulae, several choices must be made, which are listed and described in the 
following: 

1. Format of the desired image 

[0030] The format of the desired image must be chosen, i.e. the number of 
pixels in the desired image, N a , as well as the target shape of their corresponding 

SRFs, SRFj arget ( r) . Note that by specifying the target SRF shape for the pixel jt, 

one also specifies the position of the SRF's main lobe. This position will typically 
correspond to the placement of the pixel in the image space. The main lobe will then 
shift according to the pixel position. Even if the same underlying shape is chosen for 
all pixels, the target SRFs will differ by this shift. Note, however, that the target SRFs 
may also vary from pixel to pixel in terms of the underlying shape. Note also that the 
spatial arrangement of the target SRFs does not need to be regular, e.g., such that 
the main lobe positions form a Cartesian grid. Each pixel's target SRF shape may be 
chosen fully independently, e.g., such that main lobes form a hexagonal or radial 
pattern. 

[0031] In a typical embodiment, the target SRFs would be sharply localized 

functions, e.g. Dirac, boxcar, Gaussian, Lorentzian, or Sine shapes, with their 
centers arranged along a Cartesian grid. Preferably, the distance between 
neighbouring target SRF centers is chosen to be less than the smallest features of 
the encoding functions. With larger distances between SRF centers, the resolution of 
the final image will be suboptimal. 

[0032] Typically, the same target SRF shape is chosen for each pixel, with 
individual shifts corresponding to the pixel position. This has the advantage that left- 
multiplication by the matrix T, as required by Eqs. [15, 16], reduces to a convolution, 
which can be performed very efficiently using fast Fourier transform (FFT). 
[0033] Especially efficient calculations result from choosing the target SRFs to 

be Dirac shapes, with their center positions identical to the set of positions r p used 

for discretisation. In this case, the matrix Tis an identity matrix (see Eq. [9]) and can 
thus be omitted. 
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2. Discretisation 

[0034] The specific form of the discretisation must be chosen. This is done by 
specifying the total number of sampling positions, N p f and the positions 

themselves, r p . Note that these positions generally lie in 3D space. Typically, they 

should cover the object to be examined with sufficient density so that the relevant 
features of the resulting SRFs are captured in the SRF matrix according to Eq. [8]. It 
is imperative that the discretisation be sufficiently fine. In the common situation 
where Fourier encoding with gradient fields is the dominant encoding mechanism, 
the maximum distance between neighbouring sampling positions should typically not 

exceed n /(2k max ) , where k max is the maximum value that the modulus of the 
vector assumed in Eq. [2]. 

[0035] Regular Cartesian grids are a good choice for discretisation. 

Nevertheless, the positions r p may form any sampling pattern, including hexagonal 

or radial patterns, patterns with variable density, as well as patterns that are 
restricted to irregularly shaped volumes. 

3. SRF norm and spatial weighting 

[0036] In the described embodiment, SRF deviation is quantified by a vector 

norm in combination with a general linear mapping. In principle, any norm and any 
mapping may be used. A preferable choice is the Euclidean norm (or 2-norm), 
leading to the reconstruction formulae stated above. 

[0037] The additional linear mapping A can be used to introduce spatial 
weighting, thus expressing spatially varying priority in SRF optimization. Two typical 
rationales for doing so are described in the following. The first one is to exploit prior 
knowledge about the expected signal distribution in the object. Such knowledge 

permits making an estimate not just of SRF deviation but of the actual signal 
error incurred by SRF deviations. For instance, if the support of the object is known, 
A can be made diagonal, setting to zero the diagonal entries that correspond to 
positions outside the object, while the others are set to one (see J Tsao. IEEE Trans 



9 



Med Imaging 2001; 20: 1178-1183.). In this fashion, one can express that the SRF 
does not matter where no MR signal is expected and thus leave more room for noise 
minimisation. 

[0038] Often, more than just the support of the object is known. For instance, 
reference or scout images can provide an estimate of the resonance signal 
distribution inside the object. This knowledge can be exploited by setting the 
diagonals of A to the corresponding signal estimates. This would usually involve 
interpolating and potentially regridding the prior image to the grid chosen for spatial 
discretisation. The benefit is again a more accurate measure of signal error, thus 
providing specific additional room for noise minimization. 

[0039] A yet more advanced option is setting 6 = AA H equal to the expected 
signal covariance in the object. In this fashion, the spatial signal correlation present 
in the object can be exploited as well. However, this option requires significant 
statistical knowledge about the signal covariance, which will often not be available. 
[0040] As a second reason for spatial weighting, SRF shaping may require 
higher SRF fidelity in certain target regions while penalizing SRF deviations to a 
smaller degree in other regions, irrespective of prior knowledge about the specific 
signal distribution within the object. For this purpose, A should be set to a diagonal 
matrix of spatial weighting coefficients, which are high where high fidelity is desired. 

[0041] For combining prior knowledge with additional weighting, A pnor and 

A add can be chosen according to the previous considerations, respectively, yielding 
the net weighting 

A = A Qdd A prior [17] 
[0042] Finally, if neither variant of spatial weighting is desired, A should be set 

to an identity matrix. 

[0043] Apart from Euclidean and weighted Euclidean norms, the concepts of 
SRF optimization and joint SRF and noise optimization, as described above, apply 
equally to any other norm. Non-Euclidean norms will generally enhance the 
computational demands significantly. However, they can also add useful features. 
The maximum (Chebyshev) norm, for instance, judges the quality of an actual SRF 
by its maximum deviation from the ideal. This leads to so-called equi-ripple 
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deviations, where the modulus maxima have the same value along each row of the 
deviation matrix SRF discrete -T . 

4. Noise covariance matrix *P 

[0044] In order to account accurately for noise propagation into the final 

image, the noise covariance matrix should typically be constructed according to 
the paper by KP Pruessmann et al. (Magn Reson Med 2001;46:638-651), based on 
experimental noise assessment as described therein. Alternatively, W may be 
replaced by an identity matrix, thus pretending that noise in the input data is mutually 
uncorrected and of equal strength throughout. As a result, the noise estimation will 
be less accurate, shifting the balance between noise and SRF deviation towards less 
favourable results. Another alternative is to ensure actually uncorrelated noise by 
decorrelating the input data in a first step, as described in the same paper (Magn 
Reson Med 2001;46:638-651). In that case, W becomes an identity matrix and can 
be discarded without penalty. 

5. Weighting coefficients cc,fi 

[0045] The weighting coefficients a and /? , introduced in Eq. [14], express 

the relative importance of minimizing noise and SRF deviation, respectively. By 
changing their ratio, the balance between the two forms of error in the final image 
can be shifted. Practically, either factor can be held constant, while the ratio is 
adapted via the other. This remaining choice is subject to the individual imaging 
situation. Typically, the resulting noise and image artifact levels should be in the 
same order of magnitude. Preferable choices for a and /? are such that cc VP and 
ft 6 approximate the actual noise covariance in the measured data and the actual 
resonance signal covariance in the object. One problem in making these choices is 
that the strength of the resonance signal is usually not known a priori. Therefore, the 
coefficient ratio is typically determined empirically or from a preliminary 
reconstruction with non-optimal a and f$ , which provides an estimate of the signal 
strength. 
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[0046] Another option is to adjust the coefficient ratio iteratively, based on the 
actual noise and artifact levels obtained in successive reconstructions. Note that for 
any given coefficients, both error measures can be obtained from Eqs. [10], [13], 
inserting F according to Eq. [15] or [16]. In particular, with iterative adjustment of 

a, ft any desired ratio of and can be obtained on a pixel by pixel basis. 
Alternatively, a certain range of coefficient ratios may be covered by a series of 
reconstructions, providing a choice of ratios between and A^ for each pixel. A 
final image can then be obtained by selecting the most suited instance for each pixel 
value. Generally, image quality changes smoothly as a function of cc,/3 , making the 
precise choice of these parameters relatively uncritical. 

[0047] Further crucial choices independent of image reconstruction concern the 
encoding strategy. If gradient encoding is used, the pattern of sampled positions in k- 
space must be chosen. Typical choices include Cartesian, hexagonal, radial, or a 
spiral patterns. However, arbitrary /c-space patterns, including random sampling, are 
equally consistent with the present invention. 

[0048] If sensitivity encoding with one or multiple receiver coils is used, the 
choice and positioning of the coils is important. For RF encoding the choice of the 
encoding scheme is likewise essential. Any combination of these encoding 
mechanisms is consistent with the presented methods. Any further encoding 
mechanism, as long as essentially linear, can be accommodated by incorporating 
the encoding effect into the encoding functions introduced in Eq. [1]. 

Implementation 

[0049] As stated in Eq. [3], the general reconstruction formula is given by: 

p = F m [18] 
[0050] In the described embodiment, the reconstruction matrix F is calculated 
according to Eq. [15], or [16]. These two cases, labeled A and B, are distinguished 
below: 

A. If the number of rows of E is less than or equal to the number of columns of E t 
the first form of F (Eqs. [15]) is preferred as it poses the smaller inversion problem, 
leading to 
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p = TOE" (EOE H + W) + m [1 9] 

where the weighting coefficients a, (3 have been incorporated in the matrices W_ 
and 0 , respectively. Using a root decomposition (e.g. the Cholesky decomposition) 
of W - ', 

^-'=L H L t [20 ] 
[0051] Eq. [19] can be expanded and expressed equivalently as: 

p = Td(LEf ((LE)d(LE) H +iyLm ^ 

where / denotes the identity matrix. The image vector p can thus be calculated in 4 
steps: 

A1. c = Lm 

A2. c' = ((LE)6(LE) H +I) + c 

A3. c" = d(LE) H c' 
A4. p = Tc" 

[0052] The implementation of steps A1 to A4 is described below. 

B. If the number of rows of E is larger than the number of columns of £, the second 
form of F(Eq. [16]) is preferred, leading to: 

p = T(E H W- l E + d + yE H W- l m [22] 

where, again, the weighting coefficients a and have been incorporated in the 

matrices *P and 9 , respectively. Using L according to Eq. [20], this can likewise be 
expanded and expressed equivalently as: 

p = T(^LE) H (LE)+d + J {LE) H Lm [23 ] 

[0053] The image vector p can thus be calculated in 4 steps: 
B1. c = Lm 



B2. 
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B3 . c" = ((LE)"(L£)+0* )V 
B4. P = Tc " 

[0054] The implementation of steps B1 to B4 is described below. 

- Steps A1 and B1 are implemented by multiplying the corresponding data by the 

noise decorrelation matrix L. This creates a new set of data, c, whose entries 
are characterized by mutually uncorrelated noise, while each having unit noise 
level. Typically, data acquired at different points in time have no or negligible 
noise correlation. However, data acquired simultaneously with multiple receiver 
coils often exhibit noise correlation. In this case, the application of the matrix L 
reduces to a linear recombination of the data from the different receiver 
channels (see KP Pruessmann et al., Magn Reson Med 2001;46:638-651). For 
the sake of simplicity these linearly combined channels will be referred to as 
noise-decorrelated channels in the following. 

When the data are decorrelated with respect to noise, the same linear 
manipulation is performed on the encoding matrix. In the most typical case 
where noise correlation is only found among different receiver coils, this 
manipulation affects only the sensitivity-related factor in the encoding matrix. 
Multiplying the matrices L and E is then practically performed by creating linear 
recombinations of the coils' sensitivity maps. These new sensitivity maps 
correspond to virtual coils characterized by receiving mutually uncorrelated, 
unit noise. For the sake of simplicity these sensitivity maps will be referred to as 
decorrelated sensitivity maps hereafter. 

- Steps A2 and B3 involve solving an equation of the form x = A* b. In case of step 

A2, jc, A, and b represent c , ([LE)d(LE) H +/) ? and c, respectively. In 

case of step B3, Jt, A and b represent c , i^J-Ef (LE}+Q + ^ and c . 

[0055] Both steps can be solved in several ways. Here, two options are 
described. The first way is to explicitly determine the inverse A + and multiply it with 
b. The matrix inversion can be accomplished, for example, by factoring A using a 
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singular value decomposition (SVD), and by inverting the decomposed factors. The 
second way is to restate the inversion problem as A x = b and calculate a solution 
procedurally. This can be done using iterative algorithms such as the conjugate 
gradient (CG) method. The second way is the preferred embodiment, particularly for 
large matrix sizes, since it does not require explicit storage of the matrix A y but only 
requires the effect of pre-multiplication with A to be applied to intermediate vectors y. 
The latter can be done procedurally in an efficient manner by noting the following: 

- Performing the multiplication (L E)y is equivalent to replicating y by the number of 

receiver coils, multiplying each replicate with one of the decorrelated sensitivity 

maps, (sampled along r p ) t performing gridding and fast Fourier transform 

(FFT) for each virtual channel and applying regridding to determine the data at 
each sampled position of /c-space for each noise-decorrelated channel. 

- Performing the multiplication (L E) H y is equivalent to gridding /c-space data for 

each channel, applying inverse FFT for each channel, multiplying each of the 
resulting images by the complex conjugate of the corresponding noise- 
decorrelated sensitivity map, regridding along r p , and summing the results over 

the virtual channels. Further details of solving A x = b procedurally can be 
found in the paper by KP Pruessmann et al. (Magn Reson Med 2001;46:638- 
651). 

[0056] With the gridding operations, this implementation is applicable to non- 
Cartesian sampling patterns, and it additionally allows for variable density sampling, 
density weighting, and all types of /c-space shutters. 

- Steps A3 and B2 involve multiplication with (L E) H . This can be performed either 

explicitly or procedurally as described above for steps A2 and B3. Step A3 
additionally involves multiplication with 9. If 0 is a diagonal matrix, this can be 
performed efficiently as an element-wise multiplication. 

- Steps A4 and B4 involve multiplication with r. If the target voxels all have the same 

shapes and differ only by a spatial shift, this multiplication amounts to a 
convolution. Then, steps A4 and B4 can be implemented efficiently by Fourier 

transforming c " to /c-space with FFT, multiplying the transformed data by the k- 
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space representation of the target voxel shape centered in the field of view, and 
inverse Fourier transforming the result back to geometric space. If the target 
voxels are not arranged along a Cartesian grid or the spatial discretisation is 
not Cartesian, the two Fourier operations require additional gridding. 
In all situations where gridding is mentioned, the usual considerations for 
oversampling, and post-gridding apodization correction. 

Practical examples of the invention 

[0057] Images (256x256, 220 mm FOV, 3 mm thickness) were acquired from a 
phantom in a Philips Intera 1.5T scanner (Philips Medical Systems B.V., Best, 
Netherlands) using a gradient echo sequence and five surface coils placed around 
the phantom. Reconstruction was performed iteratively using the conjugate gradient 
algorithm (see Fig. 2), due to its efficiency and its ability to handle non-Cartesian 
sampling as well. Equation [3] is generally ill-conditioned, due to non-orthogonality of 
the encoding functions. In the present method the ill-conditioning has been 
overcome with Tikhonov regularization, i.e. by choosing 0 as an identity matrix. 
[0058] When the full /c-space has been acquired, aliasing is eliminated by 
Cartesian SENSE reconstruction (see Fig. 3a). However, aliasing artifacts remain 
noticeable, if /c-space is 50% zero-filled (see Fig. 3c). The artifacts are not eliminated 
by applying Cartesian SENSE at a lower resolution without zero-filling first (see Fig. 
3d), but they are eliminated completely with the reconstruction described in this 
invention, which will be referred to as the minimum-norm reconstruction (see Fig. 
3b). 

[0059] Figure 4 shows sample SRFs as obtained with minimum norm (a,c) and 
conventional (b,d) reconstruction. With the latter, the subject pixel in the lower left 
part of the phantom suffers significant contamination from an aliasing lobe just 
outside the object. This contamination is largely eliminated with the minimum norm 
approach. Contamination from the pixel's close neighbourhood is also considerably 
reduced, as evident from the line plots on the right. These plots also reveal that the 
SRF sidelobe spacing is reduced and the main lobe is minutely narrowed, 
corresponding to slight resolution improvement. These effects reflect that minimum 
norm reconstruction not only accounts for sensitivity gradients, but actually exploits 
them for favourable SRF shaping. Figure 5 shows unfiltered reconstructed images 
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along with the phantom used. With the conventional method, marked artifacts occur 
in the center, which stem from contamination from the object border. These and 
minor further artifacts are largely eliminated with minimum norm reconstruction. Note 
the visible resolution enhancement near the object border where coil sensitivity 
gradients are the steepest. 

[0060] Fig. 6 shows diagrammatically a magnetic resonance imaging system 
in which the invention is used. 

[0061] The magnetic resonance imaging system includes a set of main coils 
10 whereby a steady, uniform magnetic field is generated. The main coils are 
constructed, for example in such a manner that they enclose a tunnel-shaped 
examination space. The patient to be examined is slid on a table into this tunnel- 
shaped examination space. The magnetic resonance imaging system also includes a 
number of gradient coils 11, 12 whereby magnetic fields exhibiting spatial variations, 
notably in the form of temporary gradients in individual directions, are generated so 
as to be superposed on the uniform magnetic field. The gradient coils 11, 12 are 
connected to a controllable power supply unit 21. The gradient coils 11, 12 are 
energized by application of an electric current by means of the power supply unit 21. 
The strength, direction and duration of the gradients are controlled by control of the 
power supply unit. The magnetic resonance imaging system also includes 
transmission and receiving coils 13, 15 for generating RF excitation pulses and for 
picking up the magnetic resonance signals, respectively. The transmission coil 13 is 
preferably constructed as a body coil whereby (a part of) the object to be examined 
can be enclosed. The body coil is usually arranged in the magnetic resonance 
imaging system in such a manner that the patient 30 to be examined, being arranged 
in the magnetic resonance imaging system, is enclosed by the body coil 13. The 
body coil 13 acts as a transmission aerial for the transmission of the RF excitation 
pulses and RF refocusing pulses. Preferably, the body coil 13 involves a spatially 
uniform intensity distribution of the transmitted RF pulses. The receiving coils 15 are 
preferably surface coils 15 which are arranged on or near the body of the patient 30 
to be examined. Such surface coils 15 have a high sensitivity for the reception of 
magnetic resonance signals which is also spatially inhomogeneous. This means that 
individual surface coils 15 are mainly sensitive for magnetic resonance signals 
originating from separate directions, i.e. from separate parts in space of the body of 
the patient to be examined. The coil sensitivity profile represents the spatial 
sensitivity of the set of surface coils. The transmission coils, notably surface coils, 
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are connected to a demodulator 24 and the received magnetic resonance signals 
(MS) are demodulated by means of the demodulator 24. The demodulated magnetic 
resonance signals (DMS) are applied to a reconstruction unit. The reconstruction unit 
reconstructs the magnetic resonance image from the demodulated magnetic 
resonance signals (DMS) and on the basis of the coil sensitivity profile of the set of 
surface coils. The coil sensitivity profile has been measured in advance and is 
stored, for example electronically, in a memory unit which is included in the 
reconstruction unit. The reconstruction unit derives one or more image signals from 
the demodulated magnetic resonance signals (DMS), which image signals represent 
one or more, possibly successive magnetic resonance images. This means that the 
signal levels of the image signal of such a magnetic resonance image represent the 
brightness values of the relevant magnetic resonance image. The reconstruction unit 
25 in practice is preferably constructed as a digital image processing unit 25 which is 
programmed so as to reconstruct the magnetic resonance image from the 
demodulated magnetic resonance signals and on the basis of the coil sensitivity 
profile. The digital image processing unit 25 is notably programmed so as to execute 
the reconstruction in conformity with the present invention. The image signal from 
the reconstruction unit is applied to a monitor 26 so that the monitor can display the 
image information of the magnetic resonance image (images). It is also possible to 
store the image signal in a buffer unit 27 while awaiting further processing, for 
example printing in the form of a hard copy. 

[0062] In order to form a magnetic resonance image or a series of successive 

magnetic resonance images of the patient to be examined, the body of the patient is 
exposed to the magnetic field prevailing in the examination space. The steady, 
uniform magnetic field, i.e. the main field, orients a small excess number of the spins 
in the body of the patient to be examined in the direction of the main field. This 
generates a (small) net macroscopic magnetization in the body. These spins are, for 
example nuclear spins such as of the hydrogen nuclei (protons), but electron spins 
may also be concerned. The magnetization is locally influenced by application of the 
gradient fields. For example, the gradient coils 12 apply a selection gradient in order 
to select a more or less thin slice of the body. Subsequently, the transmission coils 
apply the RF excitation pulse to the examination space in which the part to be 
imaged of the patient to be examined is situated. The RF excitation pulse excites the 
spins in the selected slice, i.e. the net magnetization then performs a precessional 
motion about the direction of the main field. During this operation those spins are 
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excited which have a Larmor frequency within the frequency band of the RF 
excitation pulse in the main field. However, it is also very well possible to excite the 
spins in a part of the body which is much larger man such a thin slice; for example, 
the spins can be excited in a three-dimensional part which extends substantially in 
three directions in the body. After the RF excitation, the spins slowly return to their 
initial state and the macroscopic magnetization returns to its (thermal) state of 
equilibrium. The relaxing spins then emit magnetic resonance signals. Because of 
the application of a read-out gradient and a phase encoding gradient, the magnetic 
resonance signals have a plurality of frequency components which encode the 
spatial positions in, for example the selected slice. The /c-space is scanned by the 
magnetic resonance signals by application of the read-out gradients and the phase 
encoding gradients. According to the invention, the application of notably the phase 
encoding gradients results in the sub-sampling of the /c-space, relative to a 
predetermined spatial resolution of the magnetic resonance image. For example, a 
number of lines which is too small for the predetermined resolution of the magnetic 
resonance image, for example only half the number of lines, is scanned in the k- 
space. 



